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Abstract 



We numerically solve the functional differential equations (FDE's) 
of 2-particle electrodynamics, using the full electrodynamic force ob- 
tained from the retarded Lienard-Wiechert potentials and the Lorentz 
force law. In contrast, the usual formulation uses only the Coulomb 
force (scalar potential), reducing the electrodynamic 2-body problem 
to a system of ordinary differential equations (ODE's). The ODE for- 
mulation is mathematically suspect since FDE's and ODE's are known 
to be incompatible; however, the Coulomb approximation to the full 
electrodynamic force has been believed to be adequate for physics. 
We can now test this long-standing belief by comparing the FDE so- 
lution with the ODE solution, in the historically interesting case of 
the classical hydrogen atom. The solutions differ. 

A key qualitative difference is that the full force involves a 'de- 
lay' torque. Our existing code is inadequate to calculate the detailed 
interaction of the delay torque with radiative damping. However, a 
symbolic calculation provides conditions under which the delay torque 
approximately balances (3rd order) radiative damping. Thus, further 
investigations are required, and it was prematurely concluded that 
radiative damping makes the classical hydrogen atom unstable. Solu- 
tions of FDE's naturally exhibit an infinite spectrum of discrete fre- 
quencies. The conclusion is that (a) the Coulomb force is not a valid 
approximation to the full electrodynamic force, so that (b) the n-body 
interaction needs to be reformulated in various current contexts such 
as molecular dynamics. 

pacs: 03.50.De, 03.30.+p, 02.30.Ks, 02.90.+p, 87.15.Aa, 31.15.-p, 03.65.Sq, 
03.65.Ta, 95.35.+d 

keywords: many-body problem, protein dynamics, functional differ- 
ential equations, relativistic many-body problem, interpretation of quantum 
mechanics. 
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1 INTRODUCTION 



1.1 Aim 

This author had earher*^^) proposed a new model of time evolution in physics 
using mixed-type functional differential equations (FDE's), with a tilt in the 
arrow of time. This paper sets aside the notion of a 'tilt', and takes up only 
the FDE's of retarded electrodynamics. The retarded case already explicitly 
incorporates certain subtle mathematical features of electrodynamics and 
relativity noticed by Poincare, but overlooked by Einstein and subsequent 
researchers. To bring out these subtleties, this paper reports on a numeri- 
cally computed solution of FDE's of the 2-body problem of classical retarded 
electro dynamics . 

The use of the full (retarded) electrodynamic 2-particle force leads to 
the formulation of the electrodynamic 2-body problem as a system of FDE's 
that have not actually been solved earlier, numerically or otherwise, despite 
some sporadic attempts in simplified situations. '^^^ In the absence of a sys- 
tematic way to solve these FDE's, a widely used alternative has been to 
approximate the full electrodynamic 2-particle force by the Coulomb force. 
This reformulates the electrodynamic 2-body problem as an easier system 
of ODE's, which can be numerically solved with exactly the same numerical 
techniques that are used for the ODE's of the classical 2-body problem of 
Newtonian gravitation. This alternative ODE formulation of the 2-particle 
electrodynamic interaction is incorporated, for example, in models of protein 
dynamics*^^^ underlying current software such as CHARMM, wasser, amber 
etc. 

This alternative ODE formulation is, however, mathematically suspect, 
for it is known that solutions of FDE's may exhibit qualitative features im- 
possible for solutions of ODE's. On the other hand, it is believed that, from 
the viewpoint of physics, the Coulomb force is an adequate approximation 
to the full electrodynamic force. 

We can now put this long-standing belief to test: our numerical solution of 
the full-force FDE's enables us compare the two solutions in the historically 
interesting context of the classical hydrogen atom. 
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1.2 The full elect rodynamic force 

In classical electrodynamics, the force between moving charges is given by 
the Lienard-Wiechert potentials combined with the Heaviside-Lorentz force 
law. The scalar and vector (retarded) Lienard-Wiechert potentials, are given 
by the expressions*^^^ 

^(^•'*^4^(ifc^' A(M)^^^(M). (1) 
One now computes the fields E, and B by computing 

dA 

E = -Vy - — , B = V X A. (2) 



The expressions for these fields are 



(6) 



B = -R X E. (4) 

c 

When these expressions for the fields E, and B are substituted into the 
Heaviside-Lorentz force law,*^^) 

F = gi (E + vi X B ) , (5) 
the force on a charge qi moving with velocity vi is given by the expression: 



F = 



R 



47reo (R ■ u) 

H X 

c 



{ [{c^ - v'^)u + Rx {u X a)] 



R X [(c^ - v^)u + R X (u X a)]] | . 



(6) 



Here, charge qi is located at r{t) at time t, while the position of the other 
charge q at time t is given by w(t), and 



R = r{t) - w{tr). 



(7) 
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In the above expression, tr is the retarded time (the time at which the back- 
ward null cone with vertex at r(t) meets the world line w(i) of the other 
charge), and is given implicitly by the equation 

||r(t)-w(t,)|| = c{t-tr). (8) 

Further, 

u = c— — V — cR — v, (9) 

and it is understood that v = w{tr) and a = w{tr) are the velocity and 
acceleration of the charge q at the retarded time tr- A similar expression 
gives the force exerted by the charge qi on the charge q. This leads to the 
formulation of the 2-body problem as a system of FDE's. 



1.3 The Coulomb approximation 

An alternative formulation approximates the full force by the Coulomb force. 
This approximation has been justified as follows. In the Coulomb gauge, 
V ■ A = 0, the scalar potential V satisfies the Poisson equation. Neglecting 
the vector potential, the force between the two particles (i.e., the force due 
to this scalar potential) is just the Coulomb force of electrostatics, 

F = ^r, (10) 

r being the instantaneous separation of the two particles. This leads to the 
formulation of the electrodynamic 2-body problem as a system of ODE's, 
similar to the 2-particle ODE's of Newtonian gravitation. 



1.4 The difficulty: FDEs vs ODEs 

However, the vector potential cannot be so easily neglected — the possible 
justification suggested above involves an all too facile analogy between the 
scalar Lienard-Wiechert potential and electrostatics. In the language of a 
text. 

Don't be fooled, though — unlike electrostatics, V by itself doesn't 
tell you E, you need to know A as well.. . .V hj itself is not a 
physically measurable quantity — all the man in the moon can 
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measure is E, and that involves A as well. Somehow it is built 
into the vector potential, in the Coulomb gauge, that whereas 
V instantaneously reflects all changes. . . E will change only after 
sufficient time has elapsed for the 'news' to arrive. 

Briefly, the Coulomb force (based on V) acts instantaneously, while the full 
force (based on E) involves a delay (assuming retarded potentials). 

Neglecting the vector potential, hence the delay, corresponds mathemat- 
ically to approximating FDE's by ODE's. This is suspect since it is known 
that FDE's are fundamentally different from ODE's, and that solutions of 
FDE's can have qualitative features impossible for solutions of ODE's. The 
physical consequences of these differences are explained at length in the au- 
thor's previous book.*^^°) Two differences of immediate concern are the fol- 
lowing. 

(a) Given a system of ODE's of the form 

x(t) = g{t, x(t)) (11) 

a unique solution can be obtained by prescribing the initial values of x at a 
single point of time, say x(0). However, given a system of FDE 's of the form 

x(t) = fit, x(t), x(t-t,)) (12) 

and assuming these FDE's to be retarded (i.e., tr > 0), this is no longer true. 
For example, consider the FDE 

x{t)=x{t-^). (13) 

Clearly, x{t) = sin(t) and x{t) = cos(t) are solutions, and since the equation 
is linear, any linear combination x{t) = asin(t) + 6cos(t) is also a solution of 
()13|). for arbitrary constants a and b, and the values of both a and b obviously 
cannot be fixed from a knowledge of a single quantity x(0). 

In fact, from the mathematical theory of FDE's, ^^^-^ it is known that pre- 
scribing the initial values of even an infinite number of derivatives, x(0), i(0), x(0) . . . 
is not adequate. To obtain a unique solution of the retarded FDE (fT^ we 
are required to prescribe the values of x{t) over an entire interval [0, —tr], in 
the form of an initial function (j). That is, to solve retarded FDE's one needs 
to prescribe past data rather than instantaneous or 'initial' (or final) data. 

(b) Since the speed of light is so large, in many actual contexts, the delay 
tj. is likely to be extremely small {ty. < 10~^^s, for the classical hydrogen 
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1.2 


Time Asymmetry of Delay 

Collapse of solutions of a retarded equation 
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Figure 1: The figure shows three different solutions of a retarded FDE x{t) = 
h{t)x{t — \), for a suitable choice of the function b{t). Three different past histories 
prescribed for t < lead to three different solutions all of which coincide for t > 1. 
Such a phase collapse is impossible with ODE where trajectories in phase space 
can never intersect. Because of this phase collapse, FDE, unlike ODE, cannot be 
solved backward, from prescribed future data. 

atom). Under these circumstances, it is tempting to approximate a FDE of 
the form (fT^ by an ODE of the form (fTTj) by arguing that if the delay ty. is 
small, then the values of x{t — t^) can be approximately obtained from the 
values of a;(0), i;(0), x(0), 'i'(O), . . ., by means of a 'Taylor' expansion. Such a 
procedure, however, is known to be erroneous because of fundamental qual- 
itative differences between solutions of FDE's and ODE's. For example, 
unlike an ODE, which can be solved either forward or backward in time, a 
retarded FDE cannot generally be solved backward in time. Fig. 1 repro- 
duced from Ref. 1 depicts 'phase collapse' — three forward solutions of an 
FDE merge into one solution, so that a unique backward solution is impos- 
sible from data prescribed on the interval [1, 2]. Intersecting trajectories in 
phase space is a phenomenon impossible with ODE's, so that the basic clas- 
sical mechanics requirement of a phase flow breaks down.^^^^ [The origin of 
this time-asymmetry is retrospectively obvious on physical grounds, since the 
force (jni) has destroyed the underlying time symmetry of electrodynamics, by 
assuming retardation.] 

Thus, on grounds of the known mathematical differences between FDE's 
and ODE's, it is reasonable to doubt, a priori, that the full electrodynamic 
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force (jni) can be validly approximated by the Coulomb force, by neglecting 
the vector potential. 

A good way to settle this doubt is to put the matter to test by com- 
paring the Coulomb-force solution with the solutions obtained using the full 
electro dynamic force. This requires a solution of the 2-body problem with 
the full electrodynamic force, and we accordingly proceed to calculate such 
a solution. 

2 SOLVING FDE's 

2.1 Method 1 

The mathematical features of FDE's briefly recapitulated above suggest that 
retarded electrodynamics involves a 'paradigm shift' from 'classical mechan- 
ics', for we now need past data rather than instantaneous data. In a debate 
on this question at Groningen, H. D. Zeh argued against any such paradigm 
shift. Zeh maintained that there was no need for past data, since data pre- 
scribed on a spacelike hypersurface (corresponding to an instant in spacetime) 
was adequate to solve the system of Maxwell's partial differential equations 
(PDE's) together with the Heaviside-Lorentz equations of motion for each 
particle. 

Given all fields on a spacelike hypersurface, one can evolve them forward 
in time for a small region. Given all fields, in the vicinity of a single particle, 
its equations of motion reduce to a system of ODE's which can be solved in 
the usual way. We will call this method 1. We note that it implicitly involves 
the simultaneous solution of a coupled system of PDE's and ODE's. 

2.2 Method 2 

In contrast to method 1, the FDE formulation of the electrodynamic 2-body 
problem may be geometrically visualised as follows. Assuming retarded 
potentials, the force on particle 1 at a given spacetime point (a, 0) is evaluated 
as follows. One constructs the backward null cone with vertex at (a, 0), and 
determines where it intersects the world-line of particle 2, at a point {b, t), 
say (see Fig. Given the world-line of the particle 2 in a neighborhood of 
(6, t) one evaluates the resulting Lienard-Wiechert potential at (a, 0), and 
uses this retarded potential to calculate the force on particle 1 at (a, 0). This 
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(a,0), 



(b,-t) 



Figure 2: To calculate the force acting on particle 1 at (a,0), one draws the 
backward null cone with vertex at (a,0). // this intersects the world line of the 
second particle at {b,—t), one calculates the Lienard-Wiechert potential at (a, 0), 
due to the motion of particle 2 at (6, —t). Thus, determination of the present force 
on particle 1 requires a knowledge of the past motion of particle 2. 

is the force given by (jB)). The force on particle 2 at any point is calculated 
similarly. Clearly, we can solve for the motion of either particle, only if we 
are given the appropriate portions of the past world line of the other particle. 
We will call this method 2. 

2.3 Relating the two methods 

The Groningen debate brought out the following difficulty. Both the above 
methods seem to have the same underlying physical principles (Newton's 
laws of motion, possibly in a generalised form suitable for special relativity + 
Maxwell's equations). How, then, can these principles admit fundamentally 
incompatible interpretations of instantaneity and history-dependence? 

2.4 The need for past data 

This issue can be readily resolved as follows. If retarded propagators (i.e.. 
Green functions, or fundamental solutions of the wave equation, or retarded 
Lienard-Wiechert potentials) are assumed, then the field at any point x re- 
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Figure 3: Method 1 seems to requires only instantaneous data corresponding to 
fields on a Cauchy hypersurface. Assuming only two particles and retarded propa- 
gators, the fields at any point x on the hypersurface relate to past particle motions 
at Xa and Xb- As x runs over the hypersurface, the points Xa and Xh will, in gen- 
eral, cover the entire past. Thus, in the retarded case, initial data for the fields is 
the same as past data for the world lines of the particles. 

lates to particle movements in the past at the points Xa and Xf, on the re- 
spective world lines of particles a and b, where the backward null cone from 
X respectively meets the world lines of the two particles(Fig. OJ. Thus, to 
know the field at a point x we should know the particle world lines at and 
around the past points Xa, and Xb- If we do this for every point a; on a space- 
like hypersurface, the points Xa, and Xb will, in general, cover the entire past 
world lines of the two particles. Thus, in the 2-body context, given the as- 
sumption of retardation, prescribing the fields on a spacelike hypersurface is 
really equivalent to prescribing the entire past histories of the two particles. 
That is, 

instantaneous data for e.m. fields = past data for world lines of particles 

That is, the PDE+ODE method 1 only hides the underlying history-dependence 
of electrodynamics. (The above remarks need to be appropriately modified if. 
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instead of assuming retardation, one assumes advanced or mixed-type prop- 
agators. For example, if we use advanced propagators, then 'anticipation' 
should be used in place of 'history- dependence', etc.) 

2.5 FDE vs ODE+PDE 

While both methods require past data on the motion of the two particles, the 
intuitive schema underlying method 1 is currently inconvenient for the actual 
process of obtaining a solution. Formally speaking, there is no well-known 
existence and uniqueness theorem for a coupled system of PDE's-|-ODE's. 
[Separate existence theorems are, of course, known for PDE's (Maxwell's 
equations, in this case), and for ODE's (to which the particle equations 
of motion reduce, if all fields are given).] Neither is there any well known 
numerical algorithm which converts the intuitive method of iteratively solving 
coupled PDE's -I- ODE's into an actual process of calculation. ^^^^ Both formal 
proof and a numerical scheme can very likely be developed without much 
difficulty. However, neither is available as of now, so method 1 cannot, as 
of now, be used to obtain a solution of the full-force electrodynamic 2-body 
problem. 

In contrast, for FDE's there already is a formal existence and uniqueness 
theorem from past data.^^^^ Further, there are well known numerical algo- 
rithms and tested computer programs available*^^^^ (though they do not have 
all of the most desirable numerical characteristics, and the current code has 
not been formally proved to be error free). Accordingly, method 2 is currently 
the method of choice. 

Incidentally, the PDE-I-ODE method 1 seems to need data on the en- 
tire past trajectories of the two particles. This is more information than 
is usually required for the FDE method 2. [Because of the numerical stiff- 
ness of the underlying equations, in practice, one is able to numerically solve 
the classical hydrogen atom with method 2 for only short time periods, of 
the order of a femto second (10~^^ s), for which only a very short portion 
of the past history is needed.] That is, method 1, despite its appearance 
of preserving instantaneity and not needing any information from the past, 
actually seems to need more information about the past than method 2! A 
very careful analysis of the method would presumably show that the solution 
by method 1 actually uses information from only that part of the Cauchy 
hypersurface within the Cauchy horizon, so that information across the en- 
tire hypersurface is not needed, and the two methods are really equivalent. 
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At present, however, that is still a conjecture. 

3 SOLUTION OF THE ELECTRODYNAMIC 
TWO-BODY PROBLEM 

3.1 Prescription of past history and discontinuities 

How should the past history of the two particles be prescribed? Existing 
physics provides no guidelines to help answer this question: exactly like fields 
on a spacehke hypersurface in method 1, it permits us to prescribe the past 
particle motions more or less arbitrarily. 

This has two consequences worth noting. Thus, (1) mathematically, the 
past data for a FDE is not required to be a solution of the FDE, and (2) the 
mathematical theory of FDE's tells us that a discontinuity may well develop 
at the initial point where the prescribed past data joins with the solution of 
the FDE. 

From a physical point of view, this past motion may be regarded as a 
bound or constrained motion — constrained, perhaps, by additional mechani- 
cal forces — and the discontinuity at the initial point may be attributed to the 
sudden removal of the additional constraint at t = 0. These discontinuities 
propagate downstream. While it is known from the general thcory*^^^^ that 
the discontinuities of a retarded FDE arc typically smoothened out over time, 
i.e., that they move over to successively higher derivatives, we must recog- 
nize that the problem we are solving is (except in 1 dimension) technically a 
ncutral(^s) FDE ( even though only retarded Lienard-Wiechert potentials are 
being used), so that no such smoothing need take place, and the discontinu- 
ities may persist. 

The discontinuity creates a mathematical difficulty: for, in classical elec- 
trodynamics, the equations of motion are formulated assuming that the par- 
ticle trajectories are at least thrice continuously differentiable (C^). The 
numerical solution of the FDE by higher-order Runge-Kutta methods also 
assumes a similarly high level of smoothness, which assumption is not justi- 
fied if the past history is prescribed arbitrarily. 

This problem of discontinuities obviously is a suitable topic for further 
research, and in the sequel we shall set it aside and rely on the methods for 
handling discontinuities that are already incorporated into existing computer 
codes like retard and ARCHI^^^^, which handle discontinuities for example 
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by switching to low-order polynomial extrapolation or predictor-corrector 
methods near a discontinuity. Both these programs are so written that they 
permit a discontinuity even in the function (positions, velocities) at the initial 
point, though we can always arrange initial conditions so that there is a dis- 
continuity at most in the derivative (acceleration) at the initial point. We will 
prescribe past and initial data in such a way that discontinuities are confined 
to the acceleration and higher derivatives — though, like 'quantum jumps', it 
is not at present clear whether these are the only sorts of discontinuities that 
are 'physically acceptable'. The ARCHI program enables the tracking of such 
discontinuities, classified as 'hard' or 'soft' according to the theory of Wille 
and Baker. Hard discontinuities are those that propagate instantaneously 
between components, while soft discontinuities are those that propagate over 
time. In the present context, since there are no advanced interactions, all dis- 
continuities are soft, except for discontinuities between components (velocity, 
acceleration) representing derivatives. (Further, the ARCHi program has cer- 
tain refinements — such as the use of a cyclic queue to store past history — that 
are not immediately relevant.) 



3.2 The classical hydrogen atom 

To return to the physics, since our aim is to test how well the Coulomb force 
approximates the full electrodynamic force, let us consider the case of the 
classical hydrogen atom, without radiative damping. Suppose an electron 
and proton are for t < constrained to rotate rigidly in a classical circular 
2-body orbit, calculated using the Coulomb force. What will happen when 
this constraint is removed at t = 0? With the past history prescribed in this 
way, we now have a complete problem that can be solved by method 2. 



3.3 Equations of motion 



That is, we take the non-relativistic equations of motion to be given by the 
Heaviside-Lorentz force law for each particle: 



m— = g (E + V X B 

(JjL 



2 

+ 77 



3 47reoc3 df^ 



(14) 



The radiative damping term given by square brackets in (fT^ has been 
dropped, since it is not immediately relevant to our purpose of deciding 
whether the electrodynamic force (jH)) can be approximated by the Coulomb 
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force. (We know that, in the absence of radiation damping, with the Coulomb 
force, the initial 'Keplerian' orbit should remain stable.) We also take the 
mass to be constant though, relativistically, only the proper mass is constant. 
This assumption, again, is appropriate to the comparison we wish to make 
between the Coulomb force and the full electrodynamic force. 

3.4 Notation for the explicit system of equations 

A slight change of notation is helpful for the explicit calculation. We let 

K = , ^x ^ — . (15 

Aneomi 1712 

The mass ratio ^ is dimensionless, while k has dimensions of L^T~^ (=length 
in units in which c — 1). Denote the 3-dimensional trajectories of the two 
particles by Vi{t), V2{t). The delays r, f , and the vectors Ri, R2 are defined 
by the simultaneous equations 

= Rl, Ri = r2(t) -ri(t-r), Ri = ||Ri||, (16) 

^f' = Rl R2 = ri(t) -r2(t-f), i?2 = IIR2II. (17) 

Since here we consider only retarded solutions, for which r, r > we need to 
solve only the equations 

Ri-CT = 0, (18) 

R2-cf = 0, (19) 

to obtain r, f as functions of when suitable portions of the past trajectories 
ri(i), r2{t) are known. Specifically, for each given t, we compute r and f 
as the zeros of the functions Zi{t) = \\r2{t) — ri{t — t)\\ — cr, and Z2{f) = 
I |ri (t) — r2 (t — r) 1 1 — cf . The labelling has been chosen to make the equations 
below symmetrical. 

We also introduce some notation to describe the velocity and acceleration 
at the retarded times: 
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V2(t) = ^(t), V2r(t) = V2(t-f), a2r(t) = 

For the vector u used to simplify the appearance of 
expressions 



(iV2 

~dt 



(t-f). (21) 



we now have the 



ui = cRi — vir, and U2 = CR2 — ^21- (22) 
With these notations, the equations of motion take the exphcit form: 



(23) 
(24) 



(ivi 
dt 


= K 


E2 + 


Vl 

c 


dV2 

dt 


= UK 


El + 


V2 

c 



where 



E9 



El 



R2 



(R2 • U2)3 

Rl 



[U2 (C^ 
[Ui (C^ 



-'2r; 



^?r) 



+ R2 X (U2 X asr)] , (25) 
+ Rl X (ui X air)] • (26) 



(Rl • ui)3 

The shght asymmetry between the two equations ()23p and (j21I) is only 
apparent, but is helpful as follows. By convention we take the mass ratio 
H < 1, so the equation for the heavier particle is the one in which /i appears 
(and to calculate k we must use the mass of the lighter particle). For the 
purpose of computation we took fi as the mass ratio of the electron to the 
proton, and calculated k accordingly. 

From the numerical point of view, the inverse of the mass ratio is a 
measure of the numerical stiffness of the problem. Thus, the problem is 
moderately stiff, and the numerical stiffness can be reduced by solving for a 
shorter time interval. 

3.5 Actual system of equations to be solved numeri- 
cally 

The actual system of equations to be solved are the 12 retarded equations. 



m = fi{t,yi{i,y[{]{t - r)). 



(27) 
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where 



2/[i] = xi, y[2] = yi, yi3] = zi, (28) 

2/[4] = X2, y[5] = 2/2, 2/[6] = Z2. (29) 

2/[7] = vi^, 2/[8] = viy, yiQ] = vu, (30) 

2/[10] = ^^2x, 2/[ll] = V2y, y[12] = V2z- (31) 

and ri = (xi, yi, Zi), r2 = (0:2, 2/2, 2:2), Vi = (t;i^, Viy, v^) , V2 = {v2x, V2y, V2z)- 
Further, 

m = fi = 2/[i+6] for 1 < z < 6, (32) 

while for 7 < i < 12, the functions fi are to be computed using the equations 
of motion (ITll-P))). 

3.6 Initial functions 

The past history is prescribed as rigid two-body rotation around a common 
centre of mass, 

ri(t) = ri (cos cut, sin tut, 0), (33) 

r2(t) = r2(cos ut, sin tut, 0), (34) 

in a coordinate frame in which the centre of mass is fixed. From the definition 
of the centre of mass, rriiri + m2r2 = 0, so that Te/r^ = 1/fi, with = ri 
denoting the lighter particle according to the above convention for fi. Thus, 
r2 = —fJ^rQ. We took ri = tq = classical atomic radius, and u as the 
corresponding angular velocity required to balance the Coulomb force. We 
also need the past history of the velocities of the two particles, which is 
implicit in the above prescription. Explicitly, 

vi(t) = t'i(— sin cut, cos cut, 0), (35) 

V2(t) = t>2(— sin cut, cos cut, 0), (36) 
with vi = curi, t>2 = cur2 = —fiVi. 
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3.7 Choice of units for computation 

At this stage a new problem arises. The numerical values of the quantities we 
have proposed to use are as follows: fi = 5.436x10"^, k = — 2.528 x lO^m^s"^, 
ro = 5.3 X 10"^^ m, vq = a/ (/^/^o) = 2.18 x 10^ ms~\ uj = vq/tq = 
4.12 X 10^^ s^^, c = 2.99 x lO^ms"^. These numbers, spanning 27 orders of 
magnitude, present too wide a range for reliable practical computation. 

This problem can be overcome by choosing computational units appropri- 
ately. The appropriate computational units are those in which all quantities 
are represented by numbers that are as close as possible to 1. Hence, the 
base SI units, or units in which c = 1 are not the most appropriate, and 
it is preferable to take length in Angstroms or deci nano meters (1 dnm = 
]^g-io j-j^-j ^ ^T^^ time in centi femto seconds (1 cfs = 10"^'^ s), so that the nu- 
merical values of the dimensioned quantities are as follows: c ~ 30 dnm/cfs, 
n = -0.02528 (dnm)^ (cfs)-^, tq = 0.53 dnm, vq = 0.218 (dnm) {ds)-\ 
UJ = 0.412 (cfs)~^. The numbers now span only 3 orders of magnitude. This 
choice of computational units however means that the numerical solutions 
can be typically computed only for a few femto seconds. This is nevertheless 
adequate for our aims. 

3.8 The solutions 

The solutions were obtained using the retard program of Hairer et al.*-^^-* 
The results of the numerical computation are shown in the accompanying 
figures. These results establish that it is faulty to use a model based on 
the balance of instantaneous forces in a circular orbit, with the Coulomb 
force as in the classical Kepler problem. If past history is prescribed on the 
basis of this model, for a proton-electron pair, the evolution takes place as 
shown on the right of the y-axes in Fig. |3] Though the difference is not 
immediately visible, it is substantial. The differences from the classically 
expected trajectory are shown in Figs. ISl The differences are well outside 
the prescribed tolerances. 

3.9 The delay torque 

Apart from the purely numerical difference, there is a new qualitative feature. 
The full force acting on the electron has a tangential component along the 
direction of the electron velocity. Accordingly, it has a torque about the 
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Electron trajectory 




Figure 4: // the past history (t < 0) is prescribed as a Keplerian orbit, then, in 
the absence of radiation damping, there is no visible change in the electron orbit 
for the classical hydrogen atom. 

centre of mass, which torque initially accelerates the electron. With the 
Coulomb force, such a tangential component of force, or a torque, is clearly 
impossible: since the Coulomb force acts instantaneously, it always acts along 
the line through the instantaneous centre of mass. The time evolution of this 
torque is exhibited in Fig. [7| During the prescribed past history, this torque 
is a constant. Subsequently, it oscillates, so that the FDE solution oscillates 
about the ODE solution, and the entirety of the FDE orbit never strays very 
far from the entirety of the ODE orbit. 

4 HEURISTICS 

4.1 The need for heuristics 

Given the fundamental differences between FDE and ODE, it would have 
been astonishing if the above calculations had not shown up a difference from 
classical expectations. However, the above results are based on a mathemat- 
ical theory, which, despite its rigor, is not commonly known to physicists. 
Further, the results of a complex numerical calculation might well involve 
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xy-Difference 




Figure 5: The time evolution of the difference from the classical orbit of the 
electron. Time is in centi femto seconds, and distances are in deci nano meters. 
The zero difference part corresponds to the prescribed past history. The difference 
that was not visible in the previous figure, is now seen to be large. 
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Figure 6: The time evolution of the difference of only the x-coordinate of the 
electron from its classical expectation. Same scale as before. 



18 



Tangential force 



0.0005 




0.0004 


A A - 


0.0003 




0.0002 




0.0001 









0.0001 




0.0002 




0.0003 




0.0004 




0.0005 





-30 -20 -10 10 20 30 

Time 



Figure 7: The plot shows the time variation in the tangential component of the 
force exerted by the proton on the electron. This component is constant during the 
prescribed past history of the rigid rotation. Subsequently, it oscillates. Thus, the 
electron is alternately accelerated and retarded. 

some hidden error or subtle instability. Accordingly, one gains greater con- 
fidence in the numerical results if there is a way to understand the results 
intuitively and heuristically. This is not easy, given that the physicist's in- 
tuition is trained around ODE's, that are fundamentally different from the 
FDE's used in the above calculation. However, we may proceed as follows. 

4.2 The retarded inverse square force 

For an intuitive understanding of the above result, we need to approximate 
the force (jH)) by a simpler force. Any instantaneous approximation, such as 
the Coulomb force, would be mathematically unsound. Nothing we know, 
however, prevents us from approximating an FDE by a simpler FDE. That 
is, we may well approximate the full electrodynamic force by a simpler force 
while retaining the delay. 

Accordingly, let us consider a force which varies inversely as the square of 
the retarded distance. The retarded distance is the distance from the current 
position of the electron to the 'last seen' position of the proton (Fig. |H1). 
Such a force is not physically acceptable, since it cannot be relativistically 
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Figure 8: The figure shows two charged particles e and p, with masses and 
rup, rotating rigidly about a common centre of mass C , in circles of radii ri and 
r2 respectively, with an assumed uniform angular velocity uj. The classical force 
exerted hyp one acts along the line joining e to the 'last seen' position ofp, which 
differs from its instantaneous position by the delay angle 9 = uj 5t, where both 5t 
and d are uniquely fixed by the condition 6t = d/c, c being the speed of light. 



covariant; however, the simple expression for the force helps us to understand 
intuitively the consequences of a delay. 

4.3 The delay torque 

The key point here is not the exact magnitude of the force, but its direction: 
for this force acts along the line joining the current position of the electron 
to the last seen position of the proton. We can now immediately see the 
effect of even a tiny delay. Assume, to begin with, that the two particles are 
rotating rigidly in circular orbits about a common centre of mass. Because 
of the retardation, howsoever tiny, the last seen position of the proton can 
never be the same as its instantaneous position. Therefore, because of the 
retardation, howsoever tiny, the instantaneous force acting on the electron 
can not pass through the instantaneous centre of mass. The situation for 
the other charged particle is similar. These circumstances are anathema to 
classical mechanics, for it is a well known that, under such circumstances, 
angular momentum cannot be conserved. Thus, we have the astonishing 
situation that, purely under the action of internal forces, the system suffers 
a net torque. 
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The last seen position of the proton always lags behind its instantaneous 
position. Hence, for the situation of rigid rotation depicted in Fig. |S[ the 
force exerted by the proton on the electron has a positive component in the 
direction of the electron velocity. Hence, the force tends to accelerate the 
electron. 

We can use the retarded inverse- square force approximation also to heuris- 
tically try to understand the oscillation of the computed delay torque. The 
initial delay torque first accelerates both particles. The last-seen position of 
the proton starts moving towards and eventually falls on the instantaneous 
normal line to the curve being traced by the electron, so that the delay torque 
(about the instantaneous centre of rotation/curvature) is zero. The last-seen 
position of the other particle (proton) then moves 'ahead' of the instanta- 
neous normal line to the curve traced by the electron, where 'ahead' and 
'behind' is with respect to the direction in which the curve is being traced. 
The delay torque changes sign and starts decelerating the electron. If we 
neglect radiation damping, angular momentum is conserved on the average. 
This is similar to what the numerical computation with the full electrody- 
namic force showed. 

4.4 Radiation damping and stability 

There is an old argument to the effect that the introduction of radiation 
damping makes the classical hydrogen atom unstable. This argument is based 
on the Coulomb-force formulation of the electrodynamic 2-body problem. For 
the full electrodynamic force, the effect of introducing radiation reaction is 
no longer obvious, because of the existence of a delay torque which initially 
accelerates the electron. 

There are three difficulties here. First the exact form of the radiation 
damping term is not quite clear. Dirac's*^^^-* argument that the 5th and 
higher-order terms are too complicated for 'a simple thing like the electron' 
is an interesting argument, but not entirely convincing. Secondly, this author 
has pointed out that Dirac's derivation of radiation reaction approximates 
a FDE by a ODE using a 'Taylor' expansion in powers of the delay, and 
the mathematically suspect nature of this approximation is reinforced by the 
above calculation. Finally, putting in 3rd order radiation reaction may be 
expected to result in a numerically ill-conditioned problem, and appropriate 
codes for this are yet to be developed. (The DOPRI 4-5 code used implicitly 
in RETARD is not appropriate for stiff ODE's.) 
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Since the issues involved are complex, it is helpful to have an initial 
estimate about the outcome of such a rigorous calculation. For the purposes 
of such an initial estimate, about the effect of radiation damping, we again 
use the simple heuristic case of the retarded inverse square force, to determine 
whether there can be a balance of forces between the delay torque and 3rd 
order radiation damping. 

To this end, consider two particles interacting through the retarded in- 
verse square force. The force on the electron exerted by the proton is given 
by 

F = £!d. (37) 

The force acts in the direction of the 3-vector d , along which the proton is 
'last seen' by the electron. The 3-vector d may be represented by 

d = rp{t-6t) - re(t), (38) 

where Tp{t), and re(t) denote respectively the instantaneous position vectors 
of the proton and electron, respectively, at time t, and 6t is the delay, so that 
rp{t — 5t) is the 'last seen' position of the proton. 

Assuming that the two particles are in rigid rotation with constant an- 
gular velocity u, and referring back to Fig. |S1 we have, in 3-vector notation. 



Te = ri[cos utz + sin utj], (39) 

-rp{t - 5t) = r2[cos uj{t- 5t)i + sin uj{t - 5t) j] (40) 
~ r2 [cos ujt + uj5t sin ut] i 

+ r2 [sin ut — uj5t cos ujt]] , (41) 

assuming that uSt is small. This last assumption is justified since, for classical 
rigid rotation, c<j = t>o/r , while 5t ~ r/c, so that cu^t ~ Vq/c ~ 0.01. 
Hence, the delay torque on the electron is given by 

reXF = — reXd = — rir2Cj(5tk. (42) 

To a first approximation, we may use d ~ ri + r2 = (1 + e)ri, where e = 
rrLe/mp, is the ratio of the electron to proton mass. Accordingly, rir2 = erf, 
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and, further using 5t ~ d/c, the expression for the delay torque simphfies 
to 



delay ^ (1 + e)^ c ^ ^ 

On the other hand, for the 'average' force due to (3rd order) radiation 
reaction, we have the well-known formula*^^^^ 

Frad = 6^^^^ = 4^ • 3^^^ = 3^^- (^4) 
For the case under consideration, 

da 

^ = ^- --v.. (45) 
Hence, the torque due to radiation damping 



2c 2„., 2e 32- 



^rad = X F,^d = __. -c^^.r.x Ve = - (46) 

The two torques are oppositely directed. Thus, a necessary condition for 
conservation of angular momentum is 



^delay + ^rad = 0' (^7) 



which gives 



3c3 



2^ = ^3-^'^e, (48) 



or 



2 _ 3 6 3ec^ 

~ 2(l + e)2 ■ ^ ^^^^ 

Once the constraint of rigid rotation is removed, many such equilibrium 
solutions may well exist, with or without radiative damping. 

We note, incidentally, this condition for conservation of angular momen- 
tum is quite different from the classical condition for balance of forces, which 
requires = e^/rj). Both conditions can hold simultaneously only for one 
value of r which is smaller than the Bohr radius. 
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Thus, further investigations are required to determine the exact effects of 
radiative damping, and it was prematurely concluded that radiative damping 
makes the classical hydrogen atom unstable. 

4.5 Discrete spectrum and FDE's 

The discreteness of the observed hydrogen spectrum has been regarded as 
another compelling argument against classical electrodynamics. Solutions of 
FDE's can, however, admit an infinite discrete spectrum. 

Thus, for example, consider the retarded harmonic oscillator, given by 
the linear, second order, retarded FDE 

x{t) = -x{t-l), (50) 
It is easy to see that a function of the form 

x{t) = e'"', (51) 

for complex z^, is a solution of the equation (j5(J|) if and only if Zk is a solution 
of the gtiasz-polynomial equation 

z'^e^ = -1. (52) 

It is equally easy to see that the quasi-polynomial equation has 
infinitely many complex solutions Zk = Xk i: iyk, and it is known^^^^ that 
the roots are discrete, with no cluster point, and that the large magnitude 
roots are given asymptotically by the approximate expression 

Vk = 2k7T + ti{k), (53) 

Xk = -Inyk + e2(fc), (54) 

where ei(/c) 0, and €2{k) ^ as \zk\ oo. 

(Had we dropped the retardation from (j^Uj) . and applied the same proce- 
dure instead to the usual harmonic oscillator, given by the ODE x{t) = ~x{t), 
that would have led, in the well-known way, to a quadratic polynomial 
z'^ = —1 with exactly one pair of complex conjugate roots: i, and —i.) 

Since the equation (j5U|) is linear, any linear combination of these infinitely 
many oscillatory solutions of the form (|^ is again a solution. Setting aside 
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questions of convergence, it is clear that any finite linear combination of the 
form 



with Zk given by ()54|) . say, is also a solution of ()50|) . That is, in physical terms, 
the retarded harmonic oscillator, governed by the simple linear FDE ()5U|) . 
exhibits an infinite spectrum of discrete frequencies. The general solution 
is a convergent linear combination of oscillations at an infinity of discrete 
('quantized') frequencies. As in quantum mechanics, to determine a unique 
solution one needs to know an 'initial' function. If the initial function is 
prescribed arbitrarily, we can expect discontinuities. 

The above considerations remain valid for any linear FDE with constant 
coefficients and constant delay. The locally linear approximation'-^''^ sug- 
gests that such 'quantization' is also to be expected for the non-linear FDE's 
of the electrodynamic 2-body problem. [For the FDE (jl2|) this approximation 
may be obtained simply by replacing / by its first-order 'Taylor' expansion, 
and then freezing the values of the delays t,., in a neighborhood of the point 
at which we want to approximate the solution.] 

Mere discreteness of the spectrum does not, of course, mean that the 
FDE formulation, using retarded electrodynamics, will correctly reproduce 
the actual spectrum of the hydrogen atom. Equally, mere discreteness of 
the spectrum was not an adequate argument against the classical hydrogen 
atom, which was prematurely rejected on the basis of the ODE formulation 
of the 2-body problem of electrodynamics. Though a further investigation 
of FDE's may still lead to a rejection of the classical hydrogen atom for the 
right reasons, it could, on the other hand, well help to clarify the 'missing 
link' between classical and quantum mechanics, as in the structured-time 
interpretation of quantum mechanics. 



1. The Coulomb force is not a good approximation to the full electrodynamic 
force between moving charges. 

2. The classical hydrogen atom was prematurely rejected on the basis of the 
ODE formulation of the electrodynamic 2-body problem. 



n 




(55) 



k=l 



5 CONCLUSIONS 
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3. The n-body problem involving electrodynamic interactions — as formu- 
lated, for example, in current software for molecular dynamics — needs to be 
reformulated using FDE's. 
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